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Abstract 

Accurately predicting the arrival of coronal mass ejections (CMEs) at the Earth based on remote 
images is of critical significance in the study of space weather. In this paper, we make a statistical 
study of 21 Earth directed CMEs, exploring in particular the relationship between CME initial speeds 
and transit times. The initial speed of a CME is obtained by fitting the CME with the Graduated 
Cylindrical Shell model and is thus free of projection effects. We then use the drag force model to 
fit results of the transit time versus the initial speed. By adopting different drag regimes, i.e., the 
viscous, aerodynamics, and hybrid regimes, we get similar results, with the least mean estimation 
error of the hybrid model of 12.9 hours. CMEs with a propagation angle (the angle between the 
propagation direction and the Sun-Earth line) larger than its half angular width arrive at the Earth 
with an angular deviation caused by factors other than the radial solar wind drag. The drag force 
model cannot be well applied to such events. If we exclude these events in the sample, the prediction 
accuracy can be improved, i.e., the estimation error reduces to 6.8 hours. This work suggests that 
it is viable to predict the arrival time of CMEs at the Earth based on the initial parameters with a 
fairly good accuracy. Thus, it provides a method of space weather forecast of 1-5 days following the 
occurrence of CMEs. 

Subject headings: solar-terrestrial relations — Sun: coronal mass ejections 


1. INTRODUCTION 

A coronal mass ejection (CME) is a massive eruption 
of plasma threaded with magnetic field from the Sun. 
A typical CME expels 10 14-16 g of plasma and 10 29-31 
erg of kinetic energy (Howard et al. 1985) through the 
solar corona into interplanetary space. In some extreme 
cases, CMEs have speeds over 2000 km s _1 (Yashiro et 
al. 2004; Liu et al. 2014) and may reach Earth within a 
day if they are Earth directed. CMEs are known to be 
the major source of interplanetary disturbances and are 
capable of influencing geomagnetic environments with 
shock waves, ejecta, and/or magnetic clouds (cf. Gosling 
et al. 1991; Gosling 1993; Webb et al. 1994; Burlaga et 
al. 1998; Koskinen & Huttunen 2006; Zhang et al. 2007). 
Some rare events with a very strong ejecta magnetic field, 
say, Dst < —500 nT, are potential hazards to space¬ 
craft and some modern infrastructures (e.g., Carrington 
1859; Cliver & Svalgaard 2004; Siscoe et al. 2006). There¬ 
fore, understanding the CME propagation and ultimately 
making an accurate forecast of the CME arrival time are 
important issues in solar and space physics. 

Since CMEs arrive at Earth in 1-5 days after they 
launch from the Sun (e.g., Brueckner et al. 1998; Gopal- 
swamy et al. 2000) , forecasting their arrival time with an 
accuracy of a dozen of hours is possible if we well un¬ 
derstand the key factors controlling the CME propaga¬ 
tion. Usually, studying the propagation process of CMEs 
requires remote sensing images taken during the initial 
stages of CMEs and in situ observations when CMEs 
reach the Earth. However, with only the observations of 
the initial stages at hand, some simple models are usu¬ 
ally adopted to study the follow-up propagation process 
of CMEs. Gopalswamy et al. (2000) developed a con¬ 


stant acceleration or deceleration model to account for 
the findings that fast CMEs experience a deceleration 
and slow ones tend to converge to the speed of the solar 
wind. Gopalswamy et al. (2001) further improved the 
model in which the acceleration ceases before 1 AU ac¬ 
counting for the fact that slow CMEs have approximately 
the same transit time. This reduces the prediction error 
to 10.7 hours. More sophisticated approaches to CME 
kinematics consider the equation of motion governed by 
the drag force of the solar wind (Vrsnak 2001b; Vrsnak 
& Gopalswamy 2002; Borgazzi et al. 2009; Vrsnak et al. 
2010). However, as can be seen in several validation stud¬ 
ies (e.g., Owens & Cargill 2004; Colaninno et al. 2013; 
Vrsnak et al. 2014), the kinematical model yields a pre¬ 
diction error of around 10 hours based on the data set 
available, which is not significantly improved compared 
to the constant acceleration/deceleration model. Thus 
variants of the drag force model with its parameters de¬ 
termined based on different aspects have been further 
explored (Vrsnak et al. 2004; Byrne et al. 2010; Hess & 
Zhang 2014; Mostl et al. 2014). 

An important topic in space weather forecasting is the 
estimation of the speed profile of the CME during its 
propagation. In particular, for the kinematical model, 
the CME initial speed is of crucial importance. The 
speed of a CME during the initial process can either 
be measured directly from the CME fronts in coronal 
images, or be derived by fitting the images taken from 
multi-perspectives such as the twin spacecraft STEREO 
A and B. The former method gives only the speed pro¬ 
jected on the plane of sky. The latter, however, yields 
the true speed by using geometric triangulation tech¬ 
niques (Liu et al. 2010a, b) and assuming that the leading 
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edge in the images by different spacecraft corresponds to 
the same point. In the case of only one observer, some 
other geometric methods, e.g., the fixed-$ fitting (FPF) 
(Sheeley et al. 1999; Rouillard et al. 2008), the harmonic 
mean fitting (HMF) (Lugaz 2010; Mostl et al. 2011), and 
the self-similar expansion fitting (SSEF) (Davies et al. 
2012; Mostl & Davies 2013) have been adopted by fur¬ 
ther assuming a shape of the CME front. However, such 
geometric modelings consider CME propagation only in 
the ecliptic plane that may incur large errors in estimat¬ 
ing the CME initial speed. The cone model (Fisher & 
Munro 1984) and the Graduated Cylindrical Shell (GCS) 
model (Thernisien et al. 2006, 2009; Thernisien 2011) 
have also been used for a better estimation of CME 
speeds assuming a particular shape and self similarity 
of the CME. Although various models have been applied 
to case studies, previous statistical studies of a sample 
of CMEs were mostly restricted to the projection speeds 
with the drag model (e.g., Gopalswamy et al. 2000, 2001; 
Owens & Cargill 2004). It is not accurate enough, either, 
to assume a radial propagation from the source region 
since the CME propagation trajectory sometimes devi¬ 
ates from the radial direction (e.g., Vrsnak et al. 2014). 

In this paper, we present a statistical study of 21 Earth 
directed events. These events cause geomagnetic distur¬ 
bances and serve as a good sample to test the prediction 
method of CME arrival time at the Earth. The GCS 
model is applied to determine the CME initial speeds in 
a three-dimensional perspective. A fitting on the results 
of the transit time versus the initial speed with the drag 
force model shows a mean absolute error of 12.9 hours. 
In particular, 5 CMEs are identified for their angular de¬ 
flections due to mechanisms other than the solar wind 
drag, such as interactions with the background magnetic 
field or with other CMEs. The exclusion of these events 
yields a much better prediction accuracy of about 6.8 
hours. 

2. OBSERVATIONS AND DATA ANALYSIS 
2.1. The CME Sample 

CMEs have been observed for several decades but only 
in recent years is it possible to record the images of CME 
eruption simultaneously from different perspectives, such 
as those by Solar and Heliospheric Observatory (SOHO; 
Domingo et al. 1995) and Solar Terrestrial Relations Ob¬ 
servatory (STEREO; Kaiser et al. 2008). Among all the 
CMEs, the Earth directed ones are of particular interest 
since they are potential causes of geomagnetic storms. 
An Earth directed CME can be readily identified through 
the full or half halo shape (halo CME; Howard et al. 
1982) in the coronagraph of SOHO. When a CME prop¬ 
agates into the interplanetary space, it is termed as an 
ICME and can probably be observed by in situ instru¬ 
ments to record the parameters like the plasma density, 
speed, and magnetic field by spacecraft Global Geospace 
Science WIND (Acuna et al. 1995) and Advanced Com¬ 
position Explorer (ACE; Stone et al. 1998) at 1 AU when 
the CME approaches the Earth. 

The twin spacecraft STEREO-A/B are located sep¬ 
arately and away from the Sun-Earth line, with A in 
an orbit ahead of the Earth and B behind the Earth. 
On board STEREO, the Sun Earth Connection Coronal 
and Heliospheric Investigation (SECCHI; Howard et al. 


2008) contains two coronagraphs of COR1 and COR2 
that continuously take white light images in the range 
of 1.5 — 41 ?q and 2.5 — 151?0 and with a time resolu¬ 
tion of 10 and 20 minutes, respectively. Along with the 
C2 and C3 instruments of Large Angle Spectral Corona- 
graphs (LASCO, Brueckner et al. 1995) on board SOHO 
that take white light images of a field of view (FOV) up 
to 32 Rq, the three coronagraphs provide a stereoscopic 
view of the CME and can track the CME front up to 
15 Rq. 

The data from the Heliospheric Imager (HI) on board 
STEREO is also used for verification of the interplane¬ 
tary track of the CME. HI consists of two imagers. The 
FOV of HI1 is 20° x 20° centered at 14° elongation from 
the center of the Sun and HI2 has a 70° x 70° FOV cen¬ 
tered at 53°.7 from the center of the Sun. Thus, the in¬ 
struments offer the possibility to track CMEs from near 
the Sun up to 1AU. Here, we employ J map techniques 
(Sheeley et al. 1999) to track the CME propagation and 
also verify whether a CME undergoes interactions with 
other CMEs. 

In order to make a statistical study of the CME propa¬ 
gation times and test the drag force model, we employ the 
GMU CME/ICME list compiled by Phillip Hess and Jie 
Zhang (http://solar.gmu.edu/heliophysics/index. 
php/GMU_CME/ICME_List). We select the events with un¬ 
ambiguous shock fronts in the running difference images 
of COR2 and C2/C3. Our sample comprises of 21 events 
observed during 2008-2012 that are listed in Table 1. In 
the following, the 2012 October 5 event is used to demon¬ 
strate our data processing. 

2.2. GCS Model and the CME Speed 

The initial speed of a CME is determined based on the 
GCS model. First, the height of the CME apex (i.e., the 
shock front) is fitted with ray trace programs. As is shown 
below, the GCS model provides a good morphological 
approximation to the CME and a relatively robust height 
measurement. 

For a clear identification of the CME morphology, we 
use the running difference images from the three coro¬ 
nagraphs on STEREO-A/B and SOHO (Figure 1). We 
determine in total six parameters of CME including the 
longitude, latitude, tilt angle, half angular width, aspect 
ratio, and apex height that can best approximate the 
shape of the CME front. Except for the height, we as¬ 
sume that all the other five parameters do not change 
during the propagation of the CME, that is, the CME 
undergoes a self similar expansion into the space (Davies 
et al. 2012; Mostl & Davies 2013). The goodness of fit 
under a certain set of parameters is reflected by the sim¬ 
ilarity between the wire-frame fitted CME shape and the 
actual shape shown in observed images. For each event, 
we finally choose the best set of parameters that can en¬ 
sure that the goodness of fit for the time sequence of the 
event is most favorable. 

The best fit for the 2012 October 5 event (see Figure 1 
and Table 2) gives the CME height of 10.72 Rq at 05:54 
UT. Note that an automatic similarity test is conducted 
after the set of parameters is derived manually. The test 
examines the difference between the shape of the fitted 
wire frame and the observed CME profile shown in the 
images of STEREO-A/B. For this event, the similarity 
reaches a level of 68.7%. The error of the six parameters 
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TABLE 1 

List of the CME events studied in this work 


# 

onset time a 

^init b (km s 1 ) 

lon c (°) 

Iat d (°) 

<P e (°) 

ICME start f 

type® 

transit time (hr) 

1 

2008-12-12 08:37 

363 ± 23 

6 

7 

20 

2008-12-17 02:00 

EJ 

113.4 

2 

2010-04-03 09:54 

864 ±7 

6 

-25 

38 

2010-04-05 08:00 

SH±MC 

46.1 

3 

2010-04-08 03:39 

512 ± 34 

-16 

0 

28 

2010-04-11 12:00 

SH±EJ 

80.3 

4 

2010-06-16 14:39 

222 ±2 

-21 

2 

25 

2010-06-20 20:00 

EJ 

101.3 

5 

2010-12-23 05:39 

287 ±9 

18 

-22 

22 

2010-12-28 04:00 

EJ 

118.3 

6 

2011-02-15 02:24 

769 ± 12 

6 

-11 

50 

2011-02-18 03:00 

MC±CIR 

72.6 

7 

2011-03-25 05:39 

90 ±3 

-34 

0 

21 

2011-03-29 16:00 

SH±EJ 

106.3 

8 

2011-08-04 04:39 

1512 ±90 

26 

21 

67 

2011-08-05 19:00 

SH 

38.4 

9 

2011-09-06 22:54 

678 ± 13 

2 

13 

75 

2011-09-09 15:00 

SH±EJ 

64.1 

10 

2011-09-14 00:39 

505 ±5 

15 

19 

39 

2011-09-17 02:00 

SH±EJ 

73.3 

11 

2011-10-22 10:39 

882 ±4 

68 

57 

90 

2011-10-24 18:00 

SH±MC 

55.4 

12 

2011-10-27 12:39 

795 ± 52 

-36 

25 

32 

2011-11-01 08:00 

SH 

115.3 

13 

2012-01-19 14:54 

1299 ± 16 

-27 

42 

66 

2012-01-22 05:00 

SH 

62.1 

14 

2012-03-05 03:54 

1237 ± 50 

-56 

29 

63 

2012-03-08 11:00 

SH±EJ 

79.1 

15 

2012-03-13 17:39 

1616 ± 17 

51 

18 

78 

2012-03-15 13:00 

SH±EJ 

43.4 

16 

2012-03-30 15:24 

654 ±8 

-56 

30 

31 

2012-04-04 22:00 

EJ 

126.6 

17 

2012-04-19 15:39 

607 ± 15 

-24 

-30 

50 

2012-04-23 02:30 

EJ±CIR 

82.8 

18 

2012-07-12 16:39 

1224 ± 14 

0 

-11 

56 

2012-07-14 17:00 

SH±MC 

48.4 

19 

2012-09-28 00:39 

1104 ± 112 

29 

10 

57 

2012-09-30 23:00 

SH±EJ 

70.3 

20 

2012-10-05 03:39 

558 ± 21 

9 

-20 

35 

2012-10-08 05:00 

SH±MC 

73.3 

21 

2012-10-27 16:54 

340 ± 28 

12 

10 

37 

2012-10-31 15:00 

SH±MC 

94.1 


a Date and time(UT) of the first appearance in the COR2 FOV. 
b Initial speed. 

c Longitude of the propagation direction fitted by the GCS model. Earth is at 0° longitude; angles > 0° corresponds to solar west. 
d Latitude. 

e Equivalent half angular width of the CME (see Section 2.2 and 3.3 for details). 
f Arrival date and time(UT) of the associated ICME. 

s ICME properties. CIR: co-rotating interaction regions; EJ: ejecta; MC: magnetic cloud; SH: shock. 



Fig. 1.— Running difference images of the 2012 October 5 CME 
at 05:54 UT. In the bottom row, the fitted GCS model is overlaid 
as the green wire frame. 


is estimated by varying the corresponding parameter that 
results in a 10% change in similarity. The uncertainty of 
the height is of an order less than 1 Rq. Compared to the 
relatively large FOV of COR2 (15 J?q), the uncertainty 
in height is small. Thus, the height is a robust parame¬ 
ter in the GCS model fitting even though the other five 
parameters may involve relatively large uncertainties. 

A good fit yields a height-time relationship for the 
CME propagating into the space. In practice, we repeat 
the GCS fitting 6 times for each event and then make an 
average of the fitted parameters. The height-time curve 
for the 2012 October 5 event is shown in Figure 2. The 
error bars represent the standard deviation of the 6 mea¬ 
surements. The initial speed of the CME is then obtained 
through a linear fit to the curve assuming that the CME 
has already reached a static speed in the FOV of COR2 
(Zhang & Dere 2006). For this event, U;nit is 558 ±21 km 
s^ 1 . The error refers to the 3-er uncertainty of the linear 
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Fig. 2.— Height versus time curve of the 2012 October 5 event. 
The height refers to that of the CME apex obtained by the GCS 
fitting. The error bar indicates the standard deviation from 6 in¬ 
dependent measurements. The initial speed Ui n it = 558 ± 21 km 
s' 1 is obtained through a linear fit to the curve. 


fitting. 


2.3. Transit Time 

The CME transit time is the time elapsed from the 
occurrence of the CME to its arrival at 1AU. This pa¬ 
rameter can be directly obtained from the in situ data 
by WIND (see Figure 3). In the 2012 October 5 event, 
a clear shock is seen at 05:00 UT on October 8 as indi¬ 
cated by the sudden jumps in the proton density, speed, 
and temperature (as denoted by the vertical dashed line 
in the figure). Behind the shock is a sheath region with 
enhanced proton density and temperature and variable 
magnetic field. A Magnetic Cloud (MC; Burlaga et al. 
1981) is then identified by a strong magnetic field, a 
smooth rotation of the field, and a depressed proton tem¬ 
perature (marked by the shaded region). The MC ends at 
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TABLE 2 

Best fitted parameters of the 2012 October 5 CME by the GCS 

MODEL 


lon a (°) 

lat b (°) 

r c (O) 

a d (°) K e 

H f (Rq) 

54.4^ 6 0 

— 17.3+y'o 

39.1±1§ 

27.4ti» 0.410l°;“ 

in 72 - *" 0,80 

1U. o 55 


a Carrington longitude. 
b Latitude. 

c Tilt angle respect to the equator, with counterclockwise being the positive. 
d Half angular width between the two legs of the model. 
e Aspect ratio. 
f CME apex height. 
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Fig. 4.— Transit time versus initial speed for 21 CME/ICME 
events. The horizontal bars denote the 3-cr error from the linear 
fit. The four diamonds mark the events (the 2012 July 12, 2012 
October 5, 2011 October 27, and 2012 March 30 events) that serve 
as examples for CME tracking in interplanetary space with J maps. 
The 2012 October 5 event is also used as an example for the GCS 
model fitting. 

takes for the CME to propagate to the Earth, which is 
consistent with previous observations (Gopalswamy et 
al. 2001). There is also a fairly large scatter of the data 
points as shown in Figure 4. Besides the measurement 
errors in the speed and transit time, the CME may en¬ 
counter different forces (such as angular forces exerted 
by the background magnetic field and CME-CME inter¬ 
actions) as discussed below. 
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Fig. 3. — Plasma and magnetic field parameters of the ICME 
observed by WIND. The vertical dashed line shows the arrival of 
the shock while the shaded region indicates the magnetic cloud 
interval. 

17:00 UT on 2012 October 9. A connection between the 
CME in remote-sensing images and the ICME observed 
in situ can be established with the HI data in the form 
of monthly movies or J maps. Thus, the transit time is 
calculated to be the time between the first appearance of 
the CME in the COR2 FOV and the detection of shock 
near the Earth. 

3. STATISTICAL RESULTS AND PREDICTION MODELS 

Using the methods described above, we derive the ini¬ 
tial speed and the transit time for 21 CME events in 
our sample. Figure 4 plots the two parameters show¬ 
ing a close relationship between them. Note that the 
2012 October 5 event and other three events serve as 
demonstrations in the following sections. As expected, 
the faster the CME’s initial speed, the shorter time it 


3.1. The Drag Force Model 

In general, CMEs interact with the solar wind in 
three perspectives: the Lorentz force on the plasma and 
the threaded magnetic flux rope, the gravity, and the 
drag force (Chen 1996). Compared with the drag force, 
the other two can be neglected (Vrsnak & Gopalswamy 
2002), given the reduction of the magnetic field in the he¬ 
liosphere (Chen 1989; Vrsnak 1990, 2001a) and the low 
density of CMEs (Vrsnak 2001b). Thus, the drag force 
by the solar wind is a dominant force controlling the 
propagation of CMEs in the interplanetary space. 

By considering the one-dimensional problem, the equa¬ 
tion of motion can be expressed as 

^ = — j(v — w)\v — v:\P- 1 . (1) 

The exponent /3 describes the drag regime; usually /? can 
vary in the range of 1 2, with /3 = 1 for viscous and 
/3 = 2 for aerodynamics regime (Byrne et al. 2010). 
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The parameter 7 is expressed as 7 = CdAp w /m 
(Vrsnak et al. 2010), where Cd is the drag coefficient, 
m is the mass, and A is the cross section of the CME. 
We take Cd as a constant of order unity (Cargill 2004). 
Assuming that CMEs undergo a self-similar expansion 
into the space, the cross section can be calculated as 
A ~ 7r(r</?) 2 /4, where ip is the half angle, and r is 
the height of the center of the CME shell. The so¬ 
lar wind density can be given by the empirical formula 
p w ~ pi/{r/Ro) 2 (Leblanc et al. 1998), where pi is 
the empirical solar wind density at 1 AU. Therefore, 
7 ~ 7 rCdi?Q/ 9 i<^ 2 / 4 m oc p 2 /m. According to the GCS 
model, most of the CME mass is concentrated on the 
surface of the shell. Thus, p 2 /m does not change with 
CME propagation and is proportional to the initial sur¬ 
face density of the CME, which, for simplicity, we regard 
as nearly the same for different CMEs. Therefore, the 
parameter 7 is a constant in our modeling. However, in 
some previous work (e.g., Vrsnak & Gopalswamy 2002; 
Byrne et al. 2010), the parameter 7 is taken to be in¬ 
versely related to the CME height. This relation is found 
statistically from SOHO data within the FOV of LASCO 
C3 of 30 Rq Vrsnak (2001b), so it may not work for CME 
propagation in the interplanetary space. 

The equation of motion can then be integrated analyt¬ 
ically to yield 


xp 


[y(/3 - l)t + |u init - H 1_ ’] ^1, 


( 2 ) 


where 1 < f3 < 2 and the first positive sign is for Ui n i t > w 
and the first negative sign is for Ui n it < w- la particular, 
for the pure viscous (/3 = 1 ) regime, we have 

xi = wt + -(ui„it - w) (l - e -7 *) , (3) 

7 

and for the pure aerodynamic (/3 = 2 ) regime, we get 

x 2 = wt ± - In ( 7 |u in it - w\t + 1). (4) 

7 

In the follows, we refer Equation (2) to the hybrid (drag 
force) model, as discriminated from the viscous model by 
Equation (3), and the aerodynamics model by Equation 

( 4 ). 

Given a certain set of parameters 7 , /3, and w, the 
equation of motion can be solved to yield the relationship 
between the initial speed and the transit time to 1 AU. 
Thus, the theoretical transit time can be expressed as 
T = T(v i n it; 7 ,/3, w). We then get the parameters 7 , /3, 
and w by fitting the measured initial speed versus transit 
time from the sample with the theoretical relationship. 


3.2. Fitting Results for the Whole Sample 

We first fit the transit time-initial speed distribution 
for the whole sample containing the 21 CME events us¬ 
ing the drag force model, as shown in Figure 5. We im¬ 
plement the nonlinear least absolute curve fitting with 
the Levenberg-Marquardt algorithm (Markwardt 2009). 
This technique allows faster convergence to the local min¬ 
imum but dependents on the initial value of the param¬ 
eters. Considering the possible existence of several lo¬ 
cal minima in fitting errors, we choose the initial values 


from the parameter space to ensure that the final error 
is the smallest. However, after several test runs, we find 
that the mean absolute error in the prediction of transit 
time ((| Ar|)) varies little (less than 0.1 hr) subject to the 
change of the parameters. Therefore, as explained above, 
despite the mathematical difficulty, the fitted parameters 
should result in a reasonable solution. 

The parameters obtained from the fitting are listed in 
Table 3. Fitting with the hybrid model gives a mean ab¬ 
solute error (|Ar|) = 12.9 hr. Here At is defined as the 
difference between the predicted transit time and the ob¬ 
served one; therefore, At < 0 refers to cases in which the 
former is somewhat shorter than the latter. We also cal¬ 
culate the average value of the time difference, termed as 
the mean error (At), which represents a systematic de¬ 
viation of the observed transit time from the theoretical 
one. It may be caused by some physical processes that 
cannot be described in terms of the drag force model. For 
the hybrid model, (At) = —7.1 hr. The fitting by the 
pure viscous and aerodynamic models both give (| At|) = 
13.2 hr, with (At) = —5.1 hr and —9.9 hr, respectively. 

3.3. Fitting Results for a Restricted Sample 

The three-dimensional propagation of a CME is largely 
related to the ambient environment. A non-uniformly 
distributed coronal magnetic field may result in an asym¬ 
metric or non-radial motion of the CME in the early stage 
of propagation (e.g., MacQueen et al. 1986; Gopalswamy 
& Thompson 2000; Kilpua et al. 2009; Shen et al. 2011; 
Panasenco et al. 2013). A CME can also be deflected dur¬ 
ing its propagation in the interplanetary space when in¬ 
teracted with the non-radial magnetic backgrounds with 
the deflection angle dependent on the speed of the CMEs 
(e.g., Wang et al. 2004, 2014). The drag force model can¬ 
not be well applied to such CMEs with such events with 
non-radial forces. Here we denote by 9 the propagation 
angle of a CME, which is the angle between the cen¬ 
tral axis of the CME and the Sun-Earth line, and by p 
the CME half angular width. In principle, a CME with 
(p < 6 that propagates radially and self-similarly initially 
should not encounter the Earth; however, it might reach 
the Earth through angular deflections. Such events may 
worsen the fitting results and should be analyzed sepa¬ 
rately. 

The half angle p of a CME is calculated based on the 
cone model instead of the GCS model in order to re¬ 
duce the number of free parameters. Since the tilt angle 
T, half angle a , and aspect ratio n in the GCS model 
are partially degenerated and each of them involves a 
large measurement error (see Table 2 for the error derived 
from a similarity test), we combine these three parame¬ 
ters into one, which just corresponds to the half angle in 
the cone model. From the definition of the GCS model 
(Thernisien et al. 2009), we may derive that 

<P — [a(a + S)] 1/2 , (5) 

where k = sin S and T is eliminated. 

As mentioned above, those events with the criteria 
ip < 9 may undergo obvious interactions with non- 
radial forces that make their propagation direction de¬ 
viate largely from radial. After a careful examination, 
we seek out 5 events that possibly belong to this cate¬ 
gory and exclude them in the sample of CMEs. We then 
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TABLE 3 

Best fitted parameters of the drag force model 
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hybrid 

aerodynamics 
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2 
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to cn 
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546 

558+® 


<|Ar|> 

8.0 

6.8 

6.7 


(At) 

-0.7 

-3.3 

-3.5 


3 Parameters obtained by fitting the whole sample including 21 CME events. 
b Parameters obtained by fitting the restricted sample including 16 CME events. 
c The coefficient proportional to the strength of the drag force. 
d The exponent on the difference of speed of the CME and the solar wind. 
e The solar wind speed. 

f The mean absolute error between the observed transit time and the predicted one. 
s The mean error representing the systematic under- or over-estimation. 



0 500 1000 1500 2000 

Initial speed (km s _1 ) 


Fig. 5. — Fitting to the transit time-initial speed measurements for the whole sample by the drag force model. The solid black curve 
refers to the best fit by the hybrid model, with a mean absolute estimation error (|Ar|) = 12.9 hr. The dash-dotted orange curve is for the 
viscous model with (|At|) = 13.2 hr. The dashed blue curve is for the aerodynamics model with (|Ar|) = 13.2 hr. 


redo the fitting using the viscous, hybrid, and aerody¬ 
namics model to the measured results of the remainder 
16 events, as is shown in Figure 6. Note that most of the 
filtered-out events (4 out of 5) deviate largely from the 
theoretical curve, implying that the drag force models do 
not apply to them. With the restricted sample, the esti¬ 
mation errors are 8.0, 6.8, and 6.7 hours for the viscous, 
hybrid, and aerodynamics model, respectively. 

It is seen that the estimation error by the aerodynam¬ 


ics model and that by the viscous model are equally 
large when using the original sample. However, both of 
them are significantly improved if we use the restricted 
sample instead. Comparatively, the aerodynamics model 
yields an even smaller estimation error. Thus, the aero¬ 
dynamics drag might play an important role in the so¬ 
lar wind drag in the interplanetary space. The hybrid 
model, however, fails to produce a more accurate esti¬ 
mation even with an extra free parameter than the aero- 
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Fig. 6.— Fitting to the transit time-initial speed measurements for the restricted sample by the drag force model. The filled triangles 
refer to those events that are excluded from the sample (see text for details). The solid black curve refers to the best fit by the hybrid 
model, with a mean absolute estimation error (|Ar|) = 6.8 hr. The dash-dotted orange curve is for the viscous model with (|Ar|) = 9.0 
hr. The dashed blue curve is for the aerodynamics model with (|Ar|) = 6.7 hr. 


dynamics model. This is because of the existence of a 
singularity near /3 = 2 for the hybrid model (Equation 
2) where the calculated CME height diverges. Neverthe¬ 
less, the systematic error of the hybrid model is smaller 
than that of the aerodynamics model. Recently, Byrne 
et al. (2010) obtained a fitted value of the /3 parame¬ 
ter to be 2.27. Our fitting using the hybrid model gives 
a /3 value of 1.76, which seems not sufficiently large to 
support the dominant role of the aerodynamics model. 
Such a difference may be due to our use of a different 
equation of motion without the radially decreasing term. 
Judging from our statistical results and the model fitting, 
the solar wind drag can be represented by a combination 
of the viscous and aerodynamics drags, with the latter 
likely playing a major role. 

For the 5 events excluded here, the predicted tran¬ 
sit time is less than their real transit time, suggesting 
that these CMEs travel longer distances and/or at slower 
speeds than what are predicted by the model. The real 
distances are elongated by angular deflections. More¬ 
over, CME-CME interactions can slow down the CME 
propagation speed. 

3.4. Verifying the Interplanetary Propagation of CMEs 

We check in detail four representative events in our 
sample. The 2011 October 27 event has an observed 
transit time of 115.3 hours. The predicted transit time 
based on the 3 models, however, is less than 80 hours. 
Such a discrepancy suggests that the predicted propaga¬ 
tion speed is obviously higher than the actual one. The 
J map, which is constructed from the running difference 


images along the ecliptic plane by STEREO-A, is shown 
in Figure 7. The figure reveals two CMEs, indicated 
by the two arrows, the second of which is actually the 
event included in our sample, starting at 12:39 UT. One 
can obviously see that the faster one (noted as CME-2) 
catches up with the previous slower one (noted as CME- 
1 ) at a place (shown as the boxed region) where they 
undergo an interaction. Thus, the faster CME we are 
tracking should be exerted on an extra resistance force 
other than the solar wind drag. A similar process is also 
shown in the 2012 March 30 event. Another possible 
result of the interactions is the angular motion of the 
CMEs. Since in these two events the angular width of 
the CME is less than the propagation angle, interactions 
with a previous CME might cause the propagation direc¬ 
tion deviating from radial. In contrast, the 2012 July 12 
and the 2012 October 5 events do not show significant 
interplanetary interactions, so that the drag force model 
can be well applied to such cases. 

4. DISCUSSION AND CONCLUSION 

In this paper, we make statistical study of the transit 
time (between the Sun and the Earth) against the initial 
speed for a sample including 21 CME events directed 
towards the Earth. For each event, the in situ signal 
clearly shows at least one of the three identities (shock, 
ejecta, or magnetic cloud) for a typical ICME, ensuring 
that the CME does actually reach the Earth and disturb 
the geomagnetic environment. The initial speeds of these 
CMEs cover a range from 90 to 1616 km s -1 , which is 
sufficiently wide to cover both slow and fast ones. 
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Day of July 2012 


Day of October 2012 


Fig. 7.— Time-elongation maps (J maps) constructed from COR2, HU, HI2 images by STEREO A with a slice along the ecliptic plane. 
The contrast is rescaled for each row. Upper Left: The 2011 October 27 event. The two arrows point to two CME events in which CME 
2 has a higher velocity. The small box indicates the location where the two CMEs may undergo an interaction. Upper Right: The 2011 
March 30 event. The small box is the region where the two CMEs seem to merge. Lower Left: The 2012 July 12 event. Lower Right: The 
2012 October 5 event. The latter two events show no indication of interplanetary interactions. 


The initial speed of an event is determined by the GCS 
model fitting the images obtained from COR2 on board 
STEREO A/B and C2/3 on board SOHO. The images 
from three viewing angles offer a much better accuracy 
in determination of the CME height and velocity, which 
are, in particular, free from the projection effect. 

The drag force model is then applied to fit the transit 
time versus initial speed distribution. We achieve mean 
absolute errors of 13.2, 12.9, and 13.2 hours for the vis¬ 
cous, hybrid, and drag force model, respectively. These 
values are comparable with previous results (e.g., Gopal- 
swamy et al. 2001; Owens & Cargill 2004; Colaninno et al. 
2013; Vrsnak et al. 2014). However, since the drag force 
model considers only the radial drag force exerted by 
the solar wind, it cannot apply to those events that un¬ 
dergo interplanetary interactions with other CMEs and 
whose propagation direction deviates obviously from ra¬ 
dial. We pick out 5 events that likely belong to such 
a category and exclude them in the statistical sample. 
Doing so greatly improves the prediction accuracy. The 
estimation error of the transit time is reduced to be 8.0, 
6 .8, and 6.7 hours for the three models, respectively, us¬ 
ing the restricted sample. Moreover, the result shows 
that the higher order aerodynamic model achieves a bet¬ 
ter prediction than the pure viscous model. If using the 
hybrid model, the best fitted exponent is 1.76, suggest¬ 


ing that the solar wind drag is a combination of the two 
drag regimes with the aerodynamics drag somewhat in 
more effect. Since the average time for CME propaga¬ 
tion in our sample is 79.3 hours, an error of 6.7 hours 
amounts only to about 8%. Thus, it is practical to use 
the drag force model to predict the CME arrival time at 
the Earth with a fairly good accuracy. Nevertheless, to 
reach a higher accuracy, even more sophisticated mod¬ 
els are needed to account for the events with non-radial 
interactions. 

Recently, Vrsnak et al. (2014) conducted the CME ar¬ 
rival time prediction based on the drag force model fit¬ 
ting and MHD simulation. They reported estimation 
errors of about 14 hours. Here, we use the hybrid model 
and achieve a similar estimation error of 12.9 hours for 
the whole data sample but a much improved error of 6.8 
hours for the restricted sample. Besides the model, an 
accurate determination of the CME initial speeds is very 
crucial. By checking the events in our sample, we find 
that for some of them, using only the imaging observa¬ 
tions from a single instrument, assuming a radial propa¬ 
gation from the located source region, may not yield ac¬ 
curate results because of the possible angular deflections 
of the CMEs at the initial stage in the lower corona. A 
stereoscopic determination of the CME morphology and 
propagation are thus needed to deduce the CME initial 
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speed. This is the main reason why we can achieve a 
somewhat better prediction accuracy compared with pre¬ 
vious studies. Our method is based on the observations 
in the earlier process of the CME propagation, which is 
within the FOV of COR2 (15 Rq). Thus, our method 
can be used to forecast the CME arrival time once the 
event is erupted. However, a better estimation accuracy 
may be achieved by tracking the event for an extended 
distance with HI (e.g., Mostl et al. 2014). 

The three-dimensional nature of the CME morphology 
also adds to the uncertainties in the transit time estima¬ 
tion. Since CMEs have a finite angular width, those with 
a propagation direction slightly deviated from the Sun- 
Earth line may still encounter the Earth at some off-axis 
point. In this case, the CME forehead may have trav¬ 
eled a distance larger than 1 AU. Therefore, a correction 
based on the angular width and propagation angle of the 
CME is needed (Shen et al. 2014). For example, typical 


parameters of tp = 45° and 6 = 10° gives At Ri 1.8 hr, 
while ip = 60° and 6 = 30° results in At ss 15.5 hr for 
the correction to the transit time. 

Observations have also revealed that many CMEs have 
two different fronts of a bubble-shaped structure and a 
flux-rope-shaped structure (Kwon et al. 2014). In our 
study, we only model the outermost front seen in the 
running difference images. Thus, a confusion between 
the two fronts, whose speeds differ, is possible and may 
also affect the accuracy of determination of the CME 
initial speed. 


We are grateful to the referee for valuable com¬ 
ments. This work was supported by NSFC under 
grants 11303016, 11373023, and NKBRSF under grants 
2011CB811402 and 2014CB744203. 
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